#include "interpol_lin_fun.h"
#include "global.h"
#include <assert.h>

InterpolLinFun1D::InterpolLinFun1D()
{
  xlist.clear();
  ylist.clear();
  lastid = 0;
}

InterpolLinFun1D::~InterpolLinFun1D()
{
  ;
}


bool InterpolLinFun1D::load(const std::string& fname)
{
  std::ifstream inp(fname.c_str()); 
  if (!inp) {
    EOUT("Can not load file: " << fname);
    return false;
  }
  xlist.clear();
  ylist.clear();

  double x, y;
  do {
    inp >> x >> y;
    if (inp) {
      xlist.push_back(x);
      ylist.push_back(y);
    }
  } while (inp);

  inp.close();

  return true;
}

double InterpolLinFun1D::get( double x )
{
  assert(xlist.size() > 0 && ylist.size()==xlist.size());

  int    id  = lastid;
  double res = 0;

  while (1) {
    
    if (xlist[id] > x) {
      if (id == 0) { res = ylist[id]; break;} // left of first entry. constant
      id=0; 
      continue;
    }
    
    if (xlist[id] < x) {
      if (id == (int)xlist.size()-1) { res = ylist[id]; break;} // right of las entry, constant
      if (xlist[id+1] > x) {
	res = ( ( x - xlist[id] ) / ( xlist[id+1] - xlist[id] ) ) * (ylist[id+1] - ylist[id]) + ylist[id]; // linear interpolation
	break;
      }
      // xlist[id+1] <= x
      id+=1; 
      continue;
    }
    
    if (xlist[id] == x) {
      if (id == (int)xlist.size()-1) { res = ylist[id]; break;} // last entry
      if (xlist[id+1] == x) {id+=1; continue;}
      res = ylist[id];
      break;
    }
    
  }
  
  lastid = id;
  return res;
}


